rm(list = ls())

library(causalsens)
library(DirectEffects)
library(haven)
library(foreign)
library(rio)
library(tidyverse)
library(sjlabelled)
library(memisc)
library(expss)
library(dplyr)
library(tidyverse)
library(magrittr)
library(sandwich)
library(miceadds)
library(stargazer)

ian <- import("/Users/cb2257/Desktop/APSR Replication/Data/Survey/survey_final.dta")

estsample <- subset(ian, estsample == 1)

core <- c("geoid", "start_date")
estsample$treatment <- estsample$hurricane_index1 * estsample$post

estsample[core] <- lapply(estsample[core], factor)

#BLACKWELL (2014) SENSITIVITY TESTS

#CLIMATE MIGRATION ISSUE IMPORTANCE -- Figure A-8 in Main Appendix
climmig_importance_lm <- lm(climmig_importance ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climmig_importance_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

pmodel <- glm(treatment ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_pmodel <- summary(pmodel, cluster = "geoid")
summary_pmodel

alpha <- seq(-.15, .15, by = 0.01)
climmig_importance_sens <- causalsens(climmig_importance_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climmig_importance.png", width=6, height=4, units="in", res=1200)
plot(climmig_importance_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Migration Issue Importance")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climmig_importance.png", width=6, height=4, units="in", res=1200)
plot(climmig_importance_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Migration Issue Importance")
dev.off()


#CLIMATE MIGRATION POLICY ACTION -- Figure A-8 in Main Appendix
climmig_action_lm <- lm(climmig_action ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climmig_action_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

alpha <- seq(-.15, .15, by = 0.01)
climmig_action_sens <- causalsens(climmig_action_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climmig_action.png", width=6, height=4, units="in", res=1200)
plot(climmig_action_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Migration Policy Action")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climmig_action.png", width=6, height=4, units="in", res=1200)
plot(climmig_action_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Migration Policy Action")
dev.off()


#CLIMATE CHANGE ISSUE IMPORTANCE -- Figure A-8 in Main Appendix
climchg_importance_lm <- lm(climchg_importance ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_importance_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

alpha <- seq(-.15, .15, by = 0.01)
climchg_importance_sens <- causalsens(climchg_importance_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_importance.png", width=6, height=4, units="in", res=1200)
plot(climchg_importance_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Issue Importance")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_importance.png", width=6, height=4, units="in", res=1200)
plot(climchg_importance_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Issue Importance")
dev.off()


#CLIMATE CHANGE POLICY ACTION -- Figure A-8 in Main Appendix
climchg_action_lm <- lm(climchg_action ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_action_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

alpha <- seq(-.15, .15, by = 0.01)
climchg_action_sens <- causalsens(climchg_action_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_action.png", width=6, height=4, units="in", res=1200)
plot(climchg_action_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Policy Action")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_action.png", width=6, height=4, units="in", res=1200)
plot(climchg_action_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Policy Action")
dev.off()


#CLIMATE CHANGE MITIGATION POLICIES -- Figure A-8 in Main Appendix
climchg_mitigation_lm <- lm(climchg_mitigation ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_mitigation_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

alpha <- seq(-.15, .15, by = 0.01)
climchg_mitigation_sens <- causalsens(climchg_mitigation_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_mitigation.png", width=6, height=4, units="in", res=1200)
plot(climchg_mitigation_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Mitigation Policies")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_mitigation.png", width=6, height=4, units="in", res=1200)
plot(climchg_mitigation_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Mitigation Policies")
dev.off()


#CLIMATE CHANGE ADAPTATION POLICIES -- Figure A-8 in Main Appendix
climchg_adaptation_lm <- lm(climchg_adaptation ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_adaptation_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

alpha <- seq(-.15, .15, by = 0.01)
climchg_adaptation_sens <- causalsens(climchg_adaptation_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_adaptation.png", width=6, height=4, units="in", res=1200)
plot(climchg_adaptation_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Climate Change Adaptation Policies")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_adaptation.png", width=6, height=4, units="in", res=1200)
plot(climchg_adaptation_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Climate Change Adaptation Policies")
dev.off()


#SCIENCE OF CLIMATE CHANGE -- Figure A-8 in Main Appendix
climchg_science_lm <- lm(climchg_science ~ treatment + hurricane_index1 + post + republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, weights = eweight)
summary_lm <- summary(climchg_science_lm, cluster = "geoid")
covariate_summary <- summary_lm$coefficients["treatment", ]
covariate_summary

alpha <- seq(-.15, .15, by = 0.01)
climchg_science_sens <- causalsens(climchg_science_lm, pmodel, ~ republican + democrat + woman + highschool + college + age + geoid + start_date, data = estsample, alpha = alpha, confound = one.sided.att)

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sens_climchg_science.png", width=6, height=4, units="in", res=1200)
plot(climchg_science_sens,  type = "raw", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Amount of Counfounding", main="Science of Climate Change")
dev.off()

png(file="/Users/cb2257/Desktop/APSR Replication/Figures/sensr_climchg_science.png", width=6, height=4, units="in", res=1200)
plot(climchg_science_sens,  type = "r.squared", bty = "l", ylab="Effect of Hurricane Exposure", xlab="Variance Explained by Counfounding", main="Science of Climate Change")
dev.off()

